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Abstract. To determine the exact shape of the luminosity function (LF) of galaxies 
is one of the central problems in galactic astronomy and observational cosmology. 
The most popular method to estimate the LF is maximum likelihood, which is clearly 
understood with the concepts of the information theory. In the field of information 
theory and statistical inference, great advance has been made by the discovery of 
Akaike's Information Criterion (AIC). It enables us to perform a direct comparison 
among different types of models with different numbers of parameters, and provides 
us a common basis of the model adequacy. In this paper we applied AIC to the 
determination of the shape of the LF. We first treated the estimation using stepwise 
LF (Efstathiou, Ellis, & Peterson, 1988), and derived a formula to obtain the optimal 
bin number. In addition, we studied the method to compare the goodness-of-fit of 
the parametric form (Sandage, Tamman, & Yahil, 1979) with stepwise LF. 
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1. Introduction 

The luminosity function of galaxies (LF) is one of the fundamental 
descriptions of the galaxy population (e.g. Bingelli, Sandage, &: Tam- 
mann, 1988). It is also essential to interpret the galaxy number counts 
(e.g. Koo & Kron, 1992; Ellis, 1997) or to analyse galaxy clustering 
(e.g. Strauss & Willick, 1995; Efstathiou, 1996). Furthermore, the LF 
is a fundamental test for the theory of galaxy formation (e.g. Baugh, 
Cole, & Frenk, 1996). Recently, the exact shape of the LF has been of 
particular interest, because it is one of the key issues to the "faint blue 
galaxy problem" of galaxy number counts (Koo & Kron, 1992; Ellis, 
1997), and may be related to dwarf galaxy formation (e.g. Babul &; 
Rees, 1992; Babul & Ferguson, 1996 ; Hogg & Phinney, 1997; Ferguson 
& Babul, 1998). 

Instead of classical y/Vtn-estimator of the LF (Schmidt, 1968), the 
maximum likelihood method, which is free of bias induced by density 
inhomogeneity, is popular among recent studies. Various techniques 
have been proposed by expert astronomers (Lynden-Bell, 1971; Mar- 
shall et al., 1983; Choloniewski, 1986; Caditz & Petrosian, 1993). The 
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parametric method of Sandage, Tamman, & Yahil (1979) (STY) and 
the stepwise maximum likelihood method by Efstathiou, Ellis, &: Peter- 
son (1988) (EEP) are the most popular among them. And in addition, 
Lynden-Bell's method is quite sophisticated and requires no assumption 
for the probability density function. But in practice, we often need to 
smooth the obtained LF, because the discrete feature of the real data 
is not suitable for various studies. When we perform smoothing or 
binning, it remains unclear, for example, that how many numbers of 
bins we should take. Too wide bin leads to underestimation of the slope, 
and too narrow bin makes the results unstable (e.g. Caditz & Petrosian, 
1993; Strauss & Willick, 1995; Heyl et al., 1997). It is also hard to 
estimate the relative goodness between different types of statistical 
models. We should face this kind of problem, for example, when we 
try to examine whether the Schcchtcr form (Schcchtcr, 1976) provides 
an acceptable fit or not, using the likelihood ratio of the models (EEP). 
In general, the more the number of free parameter is, the better the 
fitting becomes, and consequently, the likelihood function gets larger. 

In order to evaluate the goodness-of-fit of a certain model to the 
data, the concept of the information criterion is useful. Since the middle 
of 1970's, vast advances have been made in the field of the statistical 
inference by the discovery of Akaike's Information Criterion (AIC: 
Akaike, 1974). The meaning of the AIC is clearly understood as an 
extention of the maximum likelihood method, and closely related to 
the information entropy, especially to the 'relative entropy' of two 
probability distributions. The relative entropy has a property just like 
a distance in differential geometry, i.e. it is a distance between the 
two probability distributions. Using AIC enables us to compare the 
goodness of a certain model with that of another type directly. For this 
fascinating property, AIC and its cousins are applied to the various 
fields of studies concerning statistical model selection. 

In this paper, we make an attempt to apply AIC to the estimation 
problem of the LF. In order to understand the Akaike's theory, some 
knowledges from information theory are required. Section 2 is devoted 
to this mathematical background concepts. In section 3, first wc see 
how AIC is applied to decide the step number of EEP method, and 
next how to judge the goodness of fit of the Schechter form estimated 
by STY method. Our summary is presented in section 4. 
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In this section, wc make an informal introduction of Akaike's theory. 
We do not try to be mathematically rigorous, but make an attempt to 
make it comprehensible. 



2.1. Kullback-Leibler information and information matrix 

First of all, wc consider the 'information entropy'. Here we consider the 
(discrete) probability distribution {fi}i=i,-,N- The self information is 
defined as 

/(/,) = -In/,, (1) 
where In/j = logg/j. Then, the expectation value of I{fi) is 

N 

S=E[I]=-J2fi^^fi- (2) 

1=1 

This is called the information entropy. When we have two different 
probability distributions {fi}i=i,-,N and {gi}i=i,-,N, we can construct 
the following quantity: 

^ f- 

V{f,9)^-J2fi^^-- (3) 

i=i 3' 

This is called Kullback-Leibler information (Kullback & Lcibler, 1951), 
or more comprehensively, relative entropy of {/j}t=i,-..,iv and {gi}i=i^---,N- 
For understanding the meaning of V{f, g), we define Pi = PiiOi, • " i ^k) 
{i = 1, - ■ ■ ,N), where {9i, - ■ ■ , 9k) are the parameters on which the 
model depends, and let 



< 

^gi=Pii9^^+deu---,9% + d9K) ■ 
Then, after some arithmetics, we have 

K K 

2 



(4) 



^ K K 

V{f,g) = -Y,T.^(^'^kid9kd9i, (5) 



k=l 1=1 

where 

N 



i=l 
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The subscript = 0^ means 0j = 0° for i = 1, • • • , A?^. The matrix I^i 
is called Fisher's information matrix. If wc regard I^i as a "metric" of 
the parameter space of Oj. and 0^, we can treat V{f,g) as the distance 
between the two probability distributions {/i}i=i,...,iv and {gi}i=i,-,N, 
just as we do in differential geometry. 

2.2. Maximum likelihood estimation in the context of 
information theory 

We, here, will see the maximum likelihood method considering rela- 
tive entropy V{f,g). Maximum likelihood is the method in order to 
estimate the optimal model for a given set of data. Let xi, ■ ■ ■ ,X]sr as 
a realization of random variable X which obeys the unknown prob- 
ability distribution f{X). We define a model probability distribution 
g{x;Oi, ■ ■ ■ ,9k), which depends on K parameters {dk}k=i, -,K- Then, 
the likelihood function jC is defined by 

N 

C = £{$1, ■■■ ,9k\xi,-- - ,Xn) = Y[ 9{Xn\ 0) , (7) 

n=l 

and the logarithmic likelihood is 

N 

\^C = Y,^ng{xn;0) . (8) 

n=l 

Performing maximum likelihood estimation is to find the parameter 
set 9 = 9i,---,9k which defines the most suitable model g{x;9) for 
unknown true probability distribution f{X). If we know the true f{X), 
using relative entropy eq.(3), we have 

= E[lnf{X)]-E[lngiX;9)] . (9) 

But actually we do not know f{X), we cannot obtain E [ • ]. Therefore 
we use the law of large numbers. 

THEOREM 1. (The law of large numbers) Let Xi,---,Xn inde- 
pendent random variables, which have the same mean value 

E [Xi] = • • • = E [Xm] = m , (10) 

and ^(T^ > 0, such that 

E [{Xi - rui)] <a^ for i = l,---,N . (11) 



In 



fix) 



9{X;9)_ 
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1 ^ 



i=l 



we have 



as AT — oo. 



E[X]=m, 



(12) 

(13) 



Proof of this famous theorem is seen in many textbooks on probabihty 
theory, so we omit it (see e.g. Shiryaev, 1996). 

Hence, if the size of the dataset N is sufficiently large. 



1 1 ^ 

N^''^=nY1 1^5(a:n; e) E [lng{X; 6)] 



1=1 



It leads to 



1 ^ 

V{f,g) = E[ln/(X)] - _^ln5(x„;0) 



(14) 



(15) 



i=l 



The last term can be obtained without knowing true f{X). Since 
E[ln/(X)] is independent of = 9i,---,6k, what to be maximized 
is ln£ = X^^i ln5(x„; ^) . Consequently, we obtain the likelihood 
equation: 



d_ 



\nC 



d 

ddk 



' N 

^ln5(a;„;i 



= 



for k = 1, - ■ ■ ,K, where we defined 6 to satisfy 

max(ln>C) = ln>C(^|a:;i, • • • , Xn) ■ 



(16) 



(17) 



Here, G denotes the family of parameters 6. The solution 6 of the eq. 
(16) is called the maximum likelihood estimator. 

We, next, consider the error of the logarithmic likelihood equation 
(16): 



1 ^ 

A{e) = E[ln5(X;^)]|,^^ - _Y,l^g^^^.0)\^^. . 



(18) 



n=l 



Let 6^ = 6i, ■ ■ ■ ,9^ € Q the estimator which maximize E [In ^(X; 9)]. 
Taylor expansion and theorem 1 lead the well-known important result: 



N N 



fe=i 1=1 



d9kd9i 



ln<7(X;( 



0=0° 
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N N 



(19) 



fc=i (=1 



Again, I^i is the information matrix. Thus, we can evaluate the error 
of hkelihood function using (e.g. Stuart, Ord, k, Arnold, 1999). 

2.3. Akaike's Information Criterion (AIC) 

In this subsection, we discuss Akaike's information criterion (AIC). To 
the first order of 6, likelihood equation (16) is expressed as 



= . 



(Oi - of) 

(20) 



In order to evaluate this equation, we use the multivariate central limit 
theorem. 

THEOREM 2. (Multivariate central limit theorem) Let (n = 

1, • • • , A'') a set of sample values of a vector random variable A^, which 
have a mean 



A' 



m 



(m^ < oo) and a dispersion 



(yl'=-m'=)(A'-m') 



= a 



kl 



(21) 



(22) 



for (k = 1, • • • ,K). We define as a vector random variable whose 
sample value is 



(23) 



TJien, as N —I- oo, Z'^ becomes a Gaussian vector random variable with 



Z' 



0, E, 



^k^i 



kl 



(24) 



The proof of this theorem is also incorporated in many textbooks (e.g. 
Shiryaev 1996). 

By definition of 0^, we have 



d 



lng{X-e) 



e=e° 



61=610 



= Q. (25) 
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Therefore, with eq. (6) 



e=eo 



d 



d9k 



— lng{X;e) —lng{X;e) 



d 



0=6° 



hi 



(26) 



Using the result of theorem 2, we have a vector random variable Z'^ 
which is a Gaussian with 



and whose sample value z*^ is 



n=l ^ 



e=e° 



(27) 



(28) 



We, then, apply theorem 1 to the second term of the eq.(20), it follows 
that 



If 



In g{xn;0) 



e=eo 



^2 



lng{X;e) 



6»=6»o 



= -hi- (29) 
Substituting the eqs. (28) and (29) into eq. (20) gives 

K 

z'' = VNj2^ki{ei-ef). (30) 

1=1 

If Iki is non-singular (We can always expect it by reparametrization) , 
the inverse matrix I^^^ is defined as 



K 



k=l 

where Suv is a Kronecker's delta, and it leads to 

Hence, we can write 
K K 



(31) 



(32) 



K K 



^(^) = - E il^^k - eDhm - ^?) = E E iZv^''^^ ■ (33) 



fe=i (=1 



u=l v=l 
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Now, again we use theorem 2, and eventually we obtain, 

u=l v=l 

_ _K 
~ ~N ■ 

Thus, we reach final important result: 



u=l v=l 



(34) 



1 ^ 

E[ln9iX;e)]\,^g -^In^l^nj^) 



n=l 



K 

N 



l(ln£(%i,---,x„)-i^) 



(35) 



Since, as eq. (9) shows, the distance between the true probability dis- 
tribution f{X) and model g{X; 6) is evaluated by E[ln g{X; 9)], we have 
to choose ^ G to be maximize the quantity (lni2 — K) for given K. 
Though large K provides us large ln>C, it also reduces (ln>C — K) at 
the same time, therefore the optimal number of parameters K has to 
be taken so that it maximize {lnC{6) — K). For historical reason (e.g. 
Rao, 1965), the value — 2(ln>C(^) — K) has been considered (see section 
4). It is the Akaike's information criterion (AIC). We, now, see that 
the optimal model is the one which minimizes AIC. 



3. Application of AIC to LF estimation 

3.1. The stepwise maximum likelihood method and optimal 
step number 

We apply AIC to the estimation of the shape of the LF. First we 
consider the stepwise maximum likelihood method introduced by EEP. 
Throughout this paper, for simplicity, we do not deal with some issues 
which turn out to be important when we treat the real redshift data, e.g. 
errors in magnitude measurement, surface brightness incompleteness, 
etc. These practical problems are summarized and dealt with in e.g. 
Lin et al. (1996). 

The EEP method uses the form of the LF 

K 

(l){M) = Y,<t>kW{Mk-M) (36) 

k=l 

where M is the absolute magnitude of a galaxy, which is obtained by 
M = m - 2b - b\ogdL{z) - k{z) . (37) 



taieuchi.tex; 1/02/2008; 11:21; p. 8 



9 



Here log = logio, • the apparent magnitude, dL{z) [Mpc] : the 
luminosity distance corresponding to redshift z, and k(z) is the K- 
correction. The window function W{Mi — M) is defined by 



W{Mi -M)= < 



' I Ml- AM/2 <M <Mi + AM/2 , 

(38) 

otherwise . 



In terms of information theory, the stepwise LF is regarded as a his- 
togram model, a kind of discrete probability distribution models. The 
parameters of the model 5 G 6 is in this case 9^ = (t>k {k = 1, - ■ ■ ,K) 
themselves. Given the galaxy redshift survey data of size TV, we set 

AM = ^^PPer - Mower ^ ^3g^ 

K — 1 

Mupper = max {Mi} , Miower = min {Mi} . 

i=l,---,N i=l,---,N 

The denominator of the eq. (39) is K — 1, though the number of bins 
are K, because M^ is evaluated at the center of fc-th bin, and therefore 
the magnitude range becomes Miower — AM/2 ~ Mupper + AM/2. A 
certain constraint is adopted to {(f>k}k=i,-,K in usual manner. However, 
since the normalization does not used in the determination of the shape 
of the LF (otherwise a Lagrange multiplier A would appear in the 
denominator of the eq. (44)), and neither adopted in the procedure 
of the following STY method, wc do not set any constraint in our 
formulation. According to EEP, the likelihood function is 



C{{Mk=l,-.,K\{Mi}i=i...,N) 

^ -Q E^i W{Mi - Mi)(t)i 



t=i E£i ct^iH{M^^{zi) - Mi)AM 



(40) 



H{M^in,{Zi 

1 



Ml) 



Mii^izi) - AM/2 > Ml 



Mnmizi)-Mi 1 
AM + 2 







{zi) - AM/2 < 

Ml < Mii^{zi) + AM/2 
Mii^izi) + AM/2 < Ml 



(41) 



where Mii^{zi) is the absolute magnitude corresponding to the survey 
limit miim at redshift Zi. This likelihood function C clearly depends 
on the bin width AM, and consequently, its likelihood ratio to other 
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model depends on AM. This has been regarded as an "artificial effect" 
to be eliminated by certain procedures, but it is not true because the 
choice of AM is in this case the selection of histogram model itself. 
The logarithmic likelihood is expressed as 



ln>C = 5^ ^W{Mi- Mi) In (1,1 




K 



1=1 

H{Mii^{zi) - Mi)AM 



(42) 



Hence, likelihood equation becomes 

dlnC _^W{Mk- Mi 



y H{M^^{zi) - Mk)AM _ ^ 
iriYLi4>lH{Mv,u.{zi) - Mi)AM 

and it reduces to 

A. AM EliW{Mk-Mi) 

TLi <f>iH{MiiUz,) - Mi)AM 

This equation can be solved by iteration, and we obtain the maximum 
likelihood estimator </> = {(l)k}k=i.- -,K- Thus, 

AlCEEP = -2(ln£|^^^-K) . (45) 

Therefore, the step number K should be taken so that it minimizes the 
eq. (45). 

REMARK 1. The obtained number K may not stand for the number 
of physical parameters, viz. when we get a certain K, it does not mean 
we need K physical quantities for explanation. The obtained stepwise 
LF is the one which best reflects the property of the underlying data 
population. 



3.2. Comparison of Schechter form with stepwise LF 

We, here, consider how to compare the goodness of fit of the Schechter 
form to that of EEP stepwise LF. We set the LF as 

0(M) = O.41nlO</)*(loO-4(^'-^))'^"exp(-loO-4(^'-^)) . (46) 
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The likelihood function is therefore 
/:({«, M,}|{Ma,=i,...,;v) 

TV j^iQ0.4(M,-Mi)j (^_iQOA{M,-Mi)'^ 

-1 y (^100.4(M,-M)^ (^_io0.4(M.-M)^ 

Thus the logarithmic hkelihood becomes 

N N 



InjC = 0.41nlO(l + a)^(M*-Mi)-^10°-^(^*-^') 

i=l i=l 
(•Mii^{zi) ,\ 1+a / , ,\ 

-^Iny (^io0.4(M,-M)^ + exp (-10°-4(^*-^)) dM 

= 0.41nlO(l + a)^(M* -Mi) 

i=l 

N N 

_^100-4(M.-M0_0_4ij^l0^1^r(i + „^y.^ ^ (48) 



where 



y. ^ ;^q0.4(m.-Mh^(.,)) , (49) 

/•oo 

r(z,p) = / f-^e'^dt . (50) 
Jp 

Equation (50) is known as Legcndre's incomplete gamma function. 

Since the parameters to be estimated are a and in the STY 
method, we have the set of likehhood equations: 

dlnC 



da 

= °-^^"^Q|g(^--^-)-g r(i + a,,.) 1 = °' ^''^ 

— ^ = 0.4 In 10(1 + a)iV 

-0.4 In 10 y lOO-^^* -^*) - 0.4 In 10 y " \ ^ 
ti ^r(l + a,y,)5M, 

= 0.4 In 10 
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(l + a)iV-y l00-4(A^*-^^)-0.41nl0^^ —— } 

h r(l + a,y,)/ 

= . (52) 

By solving the eqs. (51) and (52), we obtain the maximum hkehhood 
estimators a and M*. We have two free parameters in the STY method, 
and we consequently have the AIC of STY method 

AICsTY = -2(ln/:sTYle=e " 2) , (53) 

where the subscript 9 = 9 represents that a = a and = M*. 

The relative goodness of fit of the Schechter LF compared with the 
stepwise one is evaluated by 

AAIC = AICsty-AICeep 

= -2(ln/:sTY|,=^ - ln>CEEi'L=^ + K-2). (54) 



4. Summary and discussion 

In in the previous sections we introduced the Akaike's information crite- 
rion (AIC) to the maximum likehhood estimation of galaxy luminosity 
function (LF). The AIC is closely related to the "distance" between 
two probability distributions which becomes clear by using Fisher's 
information matrix. It is expressed as 

AIC = -2{\n£{9) - K) , 

where C is a likelihood function, ^ is a set of maximum likelihood 
estimators, and K is the number of free parameters of the assumed 
model. Since the concept of the information criterion seems unfamiliar 
to the astronomical community, we discuss its meaning and practical 
use in this section. 

What we must stress is that the difference between the X^-type 
goodness-of-fit quantities and AIC. The statistic usually used to eval- 
uate uncertainty of the LF estimation is the ) which was extensively 
discussed by EEP. This error evaluation is based on the well-known 
fact that the logarithmic likelinood ratio, —2\n{C{9)/C{9'^)), is asymp- 
totically distributed as distribution (e.g. Rao, 1965; Stuart, Ord, & 
Arnold, 1999). The likelihood ratio is regarded as a random variable, 
and discussed with respect to its confidence level. In order to estimate 
its distribution, EEP performed a Monte Carlo simulations and con- 
firmed their error estimation. On the other hand, as we discussed in 
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section 2, the information criterion is a value obtained as a result of 
limit theorems, and is not regarded as a random variable. In statistical 
model selection, we often use such kind of goodness-of-fit index. The 
class of information criteria including AIC was invented along with such 
concepts. Thus, though the AIC is related to the statistic, likelihood 
ratio, its value is not discussed with a confidence level (see Sakamoto, 
Ishiguro, k, Kitagawa, 1986 for details). 

Then, how should we treat the AIC value for practical use? In the 
case of stepwise LF model (EEP), if is a number of steps of the LF. 
The AIC is 

AICEEP = -2(ln/:|^^^-if), 

where = (pi,---,(pK is the step heights. We should compare the 
AIC EEP and choose the number K which minimizes AICeep- The larger 
the number of parameter is, the larger datasizc is required, because 
each parameter estimation procedure carries its own error. In case K 
is considerably larger in comparison with the data size N, we cannot 
use the limit theorems of probability theory like theorem 1, since the 
results would no longer hold. Thus, number of parameters should be 
taken as 

K<2\fN (55) 

and at most 

N 

K<-. (56) 

Bin number K comparable to N is meaninglessly large because such 
a fine binning yields horribly numerous empty bins. In such case, the 
AIC diverges by its definition, which means that the bin number is not 
a good choice. Sometimes, AIC becomes smaller and smaller as K is 
larger, and does not take minimum value. Then the assumed form of 
the model is significantly wrong, so that V{f, g) is very large. 

We next consider the comparison between the stepwise form and 
parametric models. For the STY method, free parameters are a and 
in the Schechter function (eq.(46)). In order to compare the goodness 
of fit of EEP and STY models, we can use the difference of the AIC 

AAIC = AICsTY - AICeep 

= -2(ln/:sTY|e=e - ln>CEEp|^=^ + K-2), 

where the 9 = means that a = a. and M* = M*. As we mentioned 
above, the evaluation of AIC is essentially different from that of 
statistic. Originally, EEP regarded their stepwise LF as a true proba- 
bility density function and derived the likelihood ratio. They mentioned 



takeuchi.tex; 1/02/2008; 11:21; p. 13 



14 



that the hkehhood ratio obtained from stepwise and Schechter LF was 
very large and negative. This fact is a natural consequence, because, 
from the aspect of (relative) KuUback-Leibler distance, stepwise LF 
has much larger parameters compared with parametric Schechter LF 
and therefore the goodness-of-fit is much better. Information criterion 
is suitable for such problem. 

Again we note that it is the difference of AIC values that matters 
and not the absolute values themselves. This is because we would never 
know the ture distribution from a finite size data, and we used the 
sample expectation instead of the expectation value based on true 
distribution (see eq. (15)). In our point of view, the stepwise LF is 
an estimates derived along with one of the family of statistical models 
in this paper. Thus the difference of AICs, AAIC can be used for our 
purpose. Then, what should we regard as the "scale"' of AAIC ? The 
order of the variation of the AIC is that of the number of parameters, 
K, which is an integer. Thus, when we have the difference of AICs 

AAIC>1 , (57) 

the two distributions are significantly different. Generally, AAIC = 
AICsTY — AICeep is larger than unity, therefore we judge that the 
goodness-of-fit is not sufficient, and other functional form can instead 
be used. But acutually, for the optical galaxy LF, we are interested 
in comparison between the goodness-of-fit of Schechter form and other 
forms. The functional form choice is also an interesting issue in the 
estimaiton of the IRAS galaxy LF, which is known to be significantly 
diffcrnt from Schechter form (e.g. Saunders ct al., 1996). In such case 
AIC works as a powerful tool for model selection of the fitting functions. 
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